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We generate a new complete equation of state (EOS) of nuclear matter for a wide range of 
temperatures, densities, and proton fractions ready for use in astrophysical simulations of supernovae 
and neutron star mergers. Our previous two papers tabulated the EOS at over 180,000 grid points 
in the temperature range T — to 80 MeV, the density range ub = 10 -8 to 1.6 fm -3 , and the 
proton fraction range Yp = to 0.56. In this paper we combine these data points using a suitable 
interpolation scheme to generate a single equation of state table on a finer grid. This table is 
thermodynamically consistent and conserves entropy during adiabatic compression tests. We present 
various thermodynamic quantities and the composition of matter in the new EOS, along with several 
comparisons with existing EOS tables. Our equation of state table is available for download. 

PACS numbers: 21.65.Mn,26.50.+x,26.60.Kp,21.60. Jz 



I. INTRODUCTION 



The equation of state (EOS), pressure as a function of 
density and temperature, for warm dense nuclear matter 
is an important ingredient in the theory of core collapse 
supernovae The EOS influences the shock formation 
and evolution, and determines the compactness of the 
nascent proto-neutron star. The EOS also impacts the 
position of the neutrinosphere, which is the site of last 
scattering for neutrinos in supernovae, and the emitted 
neutrino spectra [2j . This helps determine neutrino flavor 
oscillations 0], and ultimately the spectra of supernova 
neutrinos in terrestrial detectors. 

In a previous paper 0] we used a relativistic mean field 
(RMF) model to self-consistently calculate non-uniform 
matter at intermediate density and uniform matter at 
high density. In a following paper [5|, we used a Virial 
expansion for a nonideal gas of nucleons and nuclei to ob- 
tain the EOS at low densities. Our Virial approach 
uses elastic scattering phase shifts and nuclear masses as 
input. We also include Coulomb corrections that can be 
important for neutrino interactions @. Alternative ap- 
proaches to the equation of state at low densities, such 
as 10], can include model dependent strong interactions. 
Altogether our relativistic mean field and Virial EOS 
models cover the large range of temperatures, densities, 
and proton fractions shown in Table |TJ 



TABLE I: Range of temperature T, baryon density ub and 
proton fraction Yp in the EOS table. 



Parameter 


minimum 


maximum 


number of points 


T [MeV] 


0, 10" u - 8 


10 1.9 


36 


logio(ns) [fm~ 3 ] 


-8.0 


0.2 


83 


Yp 


0, 0.05 


0.56 


53 
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Note that pure neutron matter (Yp=0) is also included 
in the table. The equation of state at zero temperature 
is obtained by extrapolation from the two lowest tem- 
peratures. There are 73,840 data points from the Virial 
calculation at low densities, 17,021 data points from the 
nonuniform Hartree calculation, and 90,478 data points 
from uniform matter calculations. The overall calcula- 
tions took 7,000 CPU days in Indiana University's su- 
percomputer clusters. 

It is most efficient for astrophysical simulations of su- 
pernovae and neutron star mergers to interpolate within 
an existing EOS table to update the properties of stellar 
matter with three independent thermodynamics parame- 
ters (T, tlb and Yp). Therefore in this paper we generate, 
from these original results, a sufficiently large EOS table 
that can easily be interpolated in a way that preserves 
thermodynamic consistency. 

There exist only two realistic EOS tables that arc 
in widespread use for astrophysical simulations, the 
Lattimer-Swesty (L-S) equation of state [ll|, that uses a 
compressible liquid drop model with a Skyrme force, and 
the H. Shen, Tokh Oyamatsu and Sumiyoshi (S-S) equa- 
tion of state [TU, that uses the Thomas-Fermi and 
variational approximations with an RMF model. The 
two EOS show clear differences for example in the nu- 
clear composition during the core-collapse phase and in 
the adiabatic index around and above the phase transi- 
tion to homogeneous nuclear matter [l[. In this paper, 
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we will discuss several comparisons between our new EOS 
and these two existing ones. 

The paper is organized as follows: In section [TT| we 
present the numerical details for the interpolation scheme 
we use to generate the full EOS table. Section ITO1 shows 
results for our EOS, including various thermodynamic 
properties and the composition. The T = beta equi- 
librium EOS is also presented and used to calculate the 
neutron star mass-radius relation. We also give some 
comparisons between our EOS and the two existing EOS 
tables. Section HVl presents a summary of our results and 
gives an outlook for future work. In the appendix fVl the 
contents of our EOS table are explained in more detail 
and the public access to our EOS table is given. 

This paper presents results based on the relati vely stiff 
RMF interaction of Lalazissis, Konig, and Ring [1J| that 
has a high pressure at high densities. In later work we 
will present additional EOS tables for other softer in- 
teractions. This will allow one to more easily run astro- 
physical simulations for different EOS tables and identify 
quantities that are sensitive to the EOS. Note that a va- 
riety of uncertainties in the EOS and some laboratory 
measurements and astronomical observations related to 
the EOS have been discussed in our previous papers 0,0- 



II. NUMERICAL DETAILS 

In this section we describe how we interpolate our free 
energy results and calculate other quantities from deriva- 
tives of the free energy in a thermodynamically consistent 
manner. 



A. Equation of state at zero temperature 

At sub-MeV temperatures, the EOS is not very sensi- 
tive to variation in the temperature. In our calculations, 
the T = EOS is obtained by quadratic extrapolation 
from the results at the two lowest temperatures T = 0.158 
and 0.251 MeV |17| . Specifically, thermodynamic quan- 
tities (such as the free energy) have approximate temper- 
ature dependence, 



F(T) = F(T = 0) 



(1) 



and this equation is used to find F(T = 0). 

We start with the zero temperature T = part of the 
EOS that constitutes the most important contribution 
to the adiabatic index T and the speed of sound, par- 
ticularly at low temperature. The latter two quantities 
mostly rely on the derivative (dP/dn)\T, which involves 
a second derivative of free energy with respect to density 
and this is sensitive to noise in the bicubic interpolation 
scheme that we use. We start from the original Table HI 
which is used to generate the pressure at zero T (we also 
include the electron pressure in the following smoothing 



procedure) . The thermodynamic pressure Pth can be ob- 
tained numerically from the free energy per baryon F/A, 



P 



,( d(F/A) \ 



(2) 



T,Yv 



We smooth the pressure in a 10 points per decade ta- 
ble by removing points that differ significantly from the 
geometric mean of neighboring points and replace them 
with interpolated values. Then we use monotonic cubic 
Hermite interpolation [l5[ on the smoothed table of pres- 
sure to get a finer table with 40 points pre decade in the 
density axis. Finally we obtain the (free) energy at zero 
temperature by integrating this pressure with respect to 
density. 



B. The entropy at finite temperature 

Next, we calculate the entropy at finite temperatures 
T < 12.5 MeV. From the free energy points at tempera- 
ture T, density ns, and proton fraction Yp as in Table HI 
the entropy per baryon stu can be obtained numerically, 



Sth 



S/A 



f d(F/A) \ 



rib,Y p 



Finally, the energy per baryon eth is 



e th = F/A - Ts 



Hi ■ 



(3) 



(4) 



Note that for higher temperatures T > 12.5 MeV, mat- 
ter is uniform for all proton fractions and densities. For 
uniform matter, all thermodynamic quantities can be ob- 
tained directly from relativistic mean field calculations, 
with good thermodynamic consistency. 

The bicubic interpolation scheme we use does not guar- 
antee the monotonicity of the second derivative, e.g. 



dS_ 
df 



8 2 F 
1PT 



> 0. 



(5) 



As a result, the internal energy per baryon eth is not 
guaranteed to be monotonic either, 



de th 
dT 



> 0. 



(6) 



The bicubic interpolation with slope limiter can guaran- 
tee that the free energy is monotonic and that the first 
law of thermodynamics will be satisfied. But it is not 
sufficient to preserve monotonicities of the entropy and 
internal energy, and it can not guarantee smoothness of 
(dP/dn)\T as well. 

The following scheme is used to generate a big table 
of EOS that is thermodynamically consistent, keeping 
monotonicities of entropy and energy, and perserving the 
smoothness of (dP/dn)\T- We start with the entropy in 
the (T, rib) plane as in Table U with ten points per decade 
of density and temperature. Then we sweep through the 
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raw table of entropy, discarding numerically noisy points 
that violate the following constraints: 



dS n dS n 

ar >0 < ^ <a 



(7) 



We replace those points by interpolation from neighbor- 
ing points. Therefore we have a table of entropy that 
satisfies conditions (J7J. (The electron and photon parts 
can be handled more accurately, by monotonic interpo- 
lation on individual quantities like energy, pressure and 
entropy.) Then we further smooth the entropy as a func- 
tion of density while keeping ([7|l. by further replacing 
noisy points with interpolated values. This procedure 
ensures the smoothness of (dP/dn)\T- 

Now we perform monotonic cubic Hcrmitc interpola- 
tion on this small table of entropy with ten points per 
decade, to generate a larger EOS table with 40 points 
per decade in both the temperature and density direc- 
tions, as indicated in Table [TTJ Then we integrate this 
smoothed entropy table as a function of temperature to 
get values of the free energy (adding the energy at zero 
temperature to make the full free energy). Thus we ob- 
tain the free energy on this finer 40 by 40 points per 
decade grid in a thermodynamically consistent fashion. 



TABLE II: Range of temperature T, density ub, and proton 
fraction Yp in the finely spaced interpolated EOS table. 



Parameter 


minimum 


maximum number of points 


T [MeV] 


0, 1CT US 


10 1M 110 


logio(n s ) [fm~ 3 ] 


-8.0 


0.2 329 


Y P 


0, 0.05 


0.56 l(Y P =0)+52 



C. Bicubic Interpolation of F/A 

The final step is to carry out bicubic interpolation of 
the previous free energy values (as in Table [TTJ) to gener- 
ate the entropy and pressure by thermodynamic deriva- 
tives Eqs. (|2I3I) . This prescription guarantees the mono- 
tonicity of entropy and pressure in the final table, and 
conserves the first law of thermodynamics in adiabatic 
compression tests. 

We first apply bicubic interpolation for the free 
energy. The first derivatives on the grid points are gen- 
erated from monotonic cubic Hermite interpolation |15| . 
The second derivative - the cross derivative d 2 F/dndT, 
on the grid points is generated as in Ref . [l6| . The bicubic 



interpolation can then fit free energies with cubic func- 
tions in temperature and density coordinates, and pro- 
vides the first and second (cross) derivatives. Then the 
entropy and pressure are otained from Eqs. Il2l''>[) . Fi- 
nally, the boundary points at the highest temperature or 
density are discarded to avoid boundary artifacts in the 
interpolation. 

III. RESULTS FOR THE EQUATION OF STATE 

In this section we discuss various thermodynamic 
quantities and the composition in the EOS. First, we dis- 
cuss the free energy per baryon from Virial gas, nonuni- 
form Hartree mean field, and uniform matter calcula- 
tions. Next we map out the phase boundaries of our 
EOS. Then we show various thermodynamic quantities 
such as the pressure, entropy per baryon, chemical po- 
tentials of neutrons, protons, and electron neutrinos, and 
the adiabatic index. We also show information on the 
composition, i.e., the average mass number and proton 
number of heavy nuclei in non-uniform matter, and the 
mass fractions of different species. The zero temperature 
EOS in chemical equilibrium is presented next, and this 
is used to calculate neutron star structure. Finally, we 
make some comparisons between our EOS and those of 
Lattimer-Swesty [ll[ and H. Shen et al.fl^. 



A. Free energy and phase boundaries 

In Fig. [TJ the free energy per baryon F/A for matter 
at T = 1, 3.16, 6.31 and 10 MeV with different pro- 
ton fractions is shown as a function of density tib- The 
free energy F/A is obtained from Virial gas, nonuniform 
Hartree mean field, and uniform matter calculations, re- 
spectively. In most cases the transition (as the density 
grows) is found at the density where Hartree calcula- 
tions or uniform matter calculations give the lowest free 
energy. For matter at very low temperature (not greater 
than ~ 1 MeV) and low proton fraction (not greater than 
~ 0.1), some matching points are obtained at densities 
where Virial gas and Hartree calculations give the closest 
free energies (the difference is less than hundreds of KeV, 
and note that that mass table [l8[ used in our Virial EOS 
Q has itself 600 KeV rms deviations in nuclear binding 
energies compared to data). Detailed information has 
been reported in our previous paper |j. We reproduce 
the figure here for completeness. 



In Fig. [21 we show the phase boundaries of nuclear 
matter at different proton fractions Yp = 0.05, 0.1, 0.2 
and 0.4. The mass fraction of heavy nuclei with mass 
number A > 4 is Xa- The boundary at low densities 
indicates when Xa is greater or less than 10~ 3 . The 



boundary at high densities indicates the transition be- 
tween non-uniform matter and uniform matter. At very 
low densities, the matter is dominated by free nucleons 
and alpha particles. As the density rises, heavy nuclei 
persist to higher temperatures. Finally uniform matter 
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FIG. 1: (Color online) Free energy per nucleon of nuclear matter at temperatures of T = 1 (a), 3.16 (b), 6.31 (c), and 10 (d) 
MeV. The proton fraction ranges from Yp = 0.05 to 0.5. 



takes over at sufficiently high density. Fig. [5] also shows 
that as proton fraction rises, the temperature regime with 
appreciable heavy nuclei grows and the transition den- 
sity to uniform matter increases. The density for the 
nonuniform matter to uniform matter transition has a 
weak temperature dependence. Fig. [3] shows the phase 
diagram of nuclear matter at different temperatures T — 
1, 3.16, 6.31, and 10 MeV. Similarly as in Fig.d the left 
boundary at low density indicates where Xa is greater 



or less than 10~ 3 . The right boundary at high density 
indicates the transition between non-uniform matter and 
uniform matter. As the temperature grows, the Yp de- 
pendence of the phase boundaries becomes larger, and 
the density regime with appreciable heavy nuclei rapidly 
shrinks. Our Figs. [2] and [3] are very simillar to the cor- 
responding Figs. 2 and 3 of ref. [13| for the H Shen et 
at EOS. 



B. Pressure 

In Fig. 21 the pressure of nuclear matter with differ- 
ent proton fractions is shown as a function of density, for 
four different temperatures (1, 3.16, 6.31, and 10 MeV). 
The pressure is obtained via the bicubic Hermite inter- 
polation, in Eq. For matter at low density in the 
Virial gas phase, the pressure scales nicely with density. 



When the density increases, formation of additional nu- 
clei decreases the pressure. The Coulomb correction to 
the lattice energy reduces the pressure even further. Note 
in part of Virial gas and Hartree Wigner-Seitz cell re- 
gions, the pressure is negative. Figure 2] only shows 
hadronic contribution to the pressure. In addition one 
needs to take into account the contribution from elec- 
trons and photons, which makes the total pressure (not 



5 




FIG. 2: Phase diagram of nuclear matter at different proton fractions, (a) 0.4, (b) 0.2, (c) 0.1 and (d) 0.05. Xa is mass fraction 
of heavy nuclei with mass number A > 4. 



shown here) positive. We also note that with inclusion 
of electrons and photons in the final table, the condition 
for thermodynamic stability against density fluctuations 
is always preserved to the numerical precision of the ta- 
ble, 

dP 

— \ T > 0. (8) 
an 

Finally after the Hartree to uniform matter transition, 
the pressure rises rapidly with density because of the 
large incompressibility of nuclear mater. 

C. Entropy 

In the upper left panel of Fig. [5j the entropy per baryon 
is shown as a function of density for nuclear matter with 
T = 1 MeV at different proton fractions. The entropy is 

I 



D. Chemical potentials of n, p, and v e . 

In Fig. [6j the neutron (left subpanels, a, c, e, g) and 
proton (right subpanels, b, d, f, h) chemical potentials 
are shown for nuclear matter with Y p — 0.1 (black curves) 
and 0.5 (red curves) at different temperatures. For nu- 



calculated from Eq. ([3]). Near a density of a few times 
10~ 8 fm~ 3 , the entropy drops rapidly with density, which 
indicates formation of nuclei, first mostly alpha particles 
and then heavier nuclei at higher density. Later we will 
illustrate this by plotting the mass fractions of different 
species. For nuclear matter with a large proton frac- 
tion, larger mass fractions of heavy nuclei are present so 
that the entropy drops more rapidly with density. When 
nuclei compose more than 90% of the mass fraction, the 
entropy remains almost constant as the density increases. 
The other panels in Fig. [5] show the entropy per baryon 
in nuclear matter at higher temperatures 3.16, 6.31, and 
10 MeV. They have similar characteristics as that for T 
= 1 MeV, except that the range of densities where the 
entropy drops rapidly, due to formation of nuclei, be- 
comes smaller because it is harder to form heavy nuclei 
at higher temperatures. 



I 

clear matter at low density and high temperature, the 
chemical potentials scale with density. The formation of 
nuclei breaks the scaling behavior, i.e., the bound nucle- 
ons have much smaller chemical potentials. 

In Fig. the equilibrium electron neutrino chemical 
potential, defined as ji Ue = fi e + u p — fi n , is shown ver- 
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FIG. 3: Phase diagram of nuclear matter at different temperatures, (a) 10, (b) 6.31, (c) 3.16 and (d) 1 MeV. Xa is the same 
as Fig. |3 



sus density for nuclear matter at various temperatures 

and proton fractions. Let us focus on the top left panel 

where nuclear matter has T = 1 MeV. For symmetric 

nuclear matter, /Zj, c grows rapidly with density, since n P 

i 

almost cancels fi n and fj, Ue ~ (J-e ~ n%- For proton rich 
nuclear matter (Yp > 0.5), \i Ve will grow even faster than 
for symmetric nuclear matter. For neutron rich nuclear 
matter, /x„ e has a maximum value at some density. This 
is due to the competition between the electron chemical 
potential and the difference of the proton and neutron 
chemical potentials. The latter difference between pro- 
ton and neutron is negative when the proton fraction is 
less than 0.5. The density where /i„ c has a maximum de- 
creases with decreasing proton fraction. It is ~ 0.2 fm -3 
when Y P = 0.3, and - 0.05 fm" 3 when Y P = 0.05. In 
the limit of pure neutron matter there is no maximum 
for \i Vc , which is just negative fx n by construction. 



E. Adiabatic index T s . 

The left panel of Fig. [5] shows the temperature of the 
adiabat with entropy 5=1, for nuclear matter with pro- 
ton fraction Yp — 0.1, 0.2, 0.3, and 0.4 including baryon, 
electron, positron, and photon contributions to 5. At 
low density, matter is dominated by free nucleons in the 
gas phase. There are lots of kinetic degrees of freedom 
and one needs only a small temperature for the matter 



to have 5 = 1. When the density rises, the formation of 
heavy nuclei suppresses the entropy and it takes higher 
temperatures for the matter to reach 5=1. The en- 
tropy, in this case, comes from a fraction of free nucleons 
in a gas and the partial occupation of levels in the heavy 
nuclei. When uniform nuclear matter forms at high den- 
sity, it takes even a higher temperature to reach the same 
entropy 5=1 because the uniform matter is highly de- 
generate and the entropy comes only from excitations 
above the Fermi surfaces. 
The adiabatic index T,, 



r, = 



/<91nP\ 
\ Sinn / „ ' 



(9) 



describes the stiffness of the EOS at constant entropy. 
In the right panel of Fig. [SJ the adiabatic index is shown 
versus density for nuclear matter with constant entropy 
5 = 1 and Yp = 0.1, 0.2, 0.3, and 0.4. At subnuclear 
density, F s remains almost constant with small variation. 
It then rises rapidly at the transition from nonuniform to 
uniform matter. This characterizes the stiffening of the 
EOS due to the large nuclear incompressibility at high 
density. 



F. Average A and Z of heavy nuclei. 

Figure |H] shows the average mass number A of heavy 
nuclei (with A > 4) versus baryon density at different 
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FIG. 4: (Color online) Pressure of nuclear matter at temperatures of T = 1 (a), 3.16 (b), 6.31 (c), and 10 (d) MeV. The proton 
fraction ranges from Yp = 0.05 to 0.5. 



temperatures and proton fractions. Note that alpha par- 
ticles are not counted as heavy nuclei. Let us first look 
at the upper left panel, where T = 1 MeV. Nuclear shell 
effects give rise to several approximate plateaus in A vs 
density for each Yp, for example A = 12, ~ 50, ~ 80 
and ~ 100. Usually A is larger in matter with smaller 
Yp . There are oscillations in A in the Hartree mean field 
regime. This is due to both nuclear shell effects and small 
errors in the free energy minimization due to our using a 
finite step in the Wigner Seitz cell size The average 
mass A can be as large as 3,000 at high density before 



the final transition to uniform matter. At higher temper- 
atures, as shown in the other panels of Fig. [§[ A grows 
more rapidly to several thousand in a narrower range of 
density. As a result, the width of the plateaus in A with 
density become much shorter and finally vanish. 

The following Fig. [10] shows the average proton num- 
ber Z of heavy nuclei versus baryon density at different 
temperatures and proton fractions. The average proton 
number has very similar characteristics as for the average 
mass number, except that the differences in Z between 
different Yps are smaller. 



G. Mass fractions of nucleons and different nuclei. 

In Fig. [11] the mass fractions of free neutrons, free pro- 
tons, alpha particles and heavy nuclei are shown versus 
baryon density at four different temperatures for a proton 
fraction Y p — 0.05. The matter is extremely neutron rich, 



so the mass fraction of free neutrons is always apprecia- 
ble at any density. Note that in the Hartree mean field 
results there is only a single nucleus associated with each 
Wigner-Seitz cell. We define a nucleon level to be free 
when it has positive energy. The upper left panel is for 
T = 1 MeV. At very low densities - less than a few times 
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FIG. 5: (Color online) Entropy per nucleon of nuclear matter at temperatures of T = 1 (a), 3.16 (b), 6.31 (c), and 10 (d) MeV. 
The proton fraction ranges from Yp = 0.05 to 0.5. 



1CP 8 fm~ 3 , the system is basically a free neutron and pro- 
ton gas. The alpha fraction rises rapidly to 10% around a 
density of 10~ 7 fm -3 . Some of the free neutrons are used 
to form alphas, so the fraction of free neutrons drops to 
90%. The free neutron fraction drops further beyond a 
density of 10 -6 fm -3 , where heavy nuclei appear with a 
mass fraction of about 15%. The fraction of heavy nuclei 
grows faster after the Virial gas - Hartree mean field tran- 
sition around 3T0~ 4 fm~ 3 . More and more free neutrons 
are bound to form heavier nuclei, as can be seen in Fig. [5] 
The fraction of heavy nuclei reaches a maximum value of 
70% before the final transition to uniform matter. The 
other panels give the mass fractions of different species 
at higher temperatures. In these cases, the appearance 
of alpha particles and heavy nuclei occur at higher den- 
sities and with smaller mass fractions compared to the 
T = 1 MeV case. When T = 6.31 MeV, there is only 
a narrow range of densities where alpha particles have a 
mass fraction of a couple of percent and the amount of 
heavy nuclei is very small. 

Figs . IT21 and [T31 show the mass fractions as were shown 



in Fig. 111! but for matter with higher proton fractions 
Y p = 0.2 and 0.5, respectively. Alpha particles appear at 
similar densities as for Y p = 0.05, but with much larger 
mass fractions. For example, the alpha fraction can al- 
most reach 100% around 10 -6 fm~ 3 when T = 1 MeV 
and Yp = 0.5. Heavy nuclei also appear at similar densi- 
ties as in Y p = 0.05, and with much larger mass fractions. 
As the proton fraction rises, heavy nuclei can survive to 
higher temperatures, even though only at very high den- 
sity. For example, when T — 10 MeV and Yp = 0.5, heavy 
nuclei have a large mass fraction for densities around a 
few times 10 -2 fm -3 . 



H. Beta equilibrium at zero temperature. 

Matter in beta equilibrium at zero temperature is par- 
ticularly interesting and determines neutron star struc- 
ture, including the mass-radius relation and the neutron 
star crust. The T = EOS is included in our table [TTl 
which can be used to find the beta equilibrium EOS by 
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FIG. 6: (Color online) Chemical potentials of neutron and proton at different temperature and proton fractions. 



imposing the constraint fi n — fi p + fi e . The resulting en- 
ergy per baryon, pressure, and proton fraction are shown 
as function of density in Fig. I14I The energy per baryon 
starts near -8.5 MeV, a typical value for stable nuclei, 
at 1CP 8 fm~ 3 , becomes positive around 10~ 3 fm -3 , and 
finally rises rapidly after the transition to bulk uniform 
matter around 0.05 fm -3 . The pressure vs density is 
polytropic like below 10~ 4 fm -3 , and then it develops 
small plateaus. After the transition to uniform matter, 
the pressure rises more rapidly which indicates a stiff- 
ening of the EOS. At 10 ~ 8 fm~ 3 matter is made up of 
stable nuclei with similar numbers of neutrons and pro- 
tons. That's why the proton fraction at 10 -8 fm~ 3 is as 
high as 0.45. The proton fraction drops as the density 
increases, which means matter becomes more and more 
neutron rich due to the large electron Fermi energy. The 
proton fraction reaches a minimum value of 0.007 around 
0.025 fm -3 . However, it increases again due to the large 
symmetry energy at high density for the NL3 |14| interac- 
tion, which favors more nearly equal neutron and proton 
fractions. 



I. Neutron star structure 



The mass and radius of a neutron star are obtained by 
solving the Tolman-Oppenheimer-Volkoff (TOV) equa- 
tions pi, [11 . 



By using the EOS in beta equilibrium at zero tem- 
perature from the previous section in the core and using 
Baym-Pethick-Sutherland EOS in the crust [2l[, we ob- 
tain the neutron star mass-radius and mass-central den- 
sity relations, which are plotted in Fig. [15j The neutron 
star's maximum gravitational mass is about 2.77 solar 
mass with a radius of 13.3 km. The corresponding cen- 
tral density is about 10 15 g/cm 3 . The large value for the 
maximum neutron star mass is not too surprising since 
NL3 provides a very stiff EOS as we have emphasized 
several times. 
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FIG. 7: (Color online) Equilibrium electron neutrino chemical potential /i„ e for nuclear matter at various temperatures and 
proton fractions. It is defined as /i„ e = fi e + fi p — fi„. 




FIG. 8: (Color online) Left: temperature along adiabat S = 1, for nuclear matter, including baryon, electron, positron and 
photon, and proton fraction Yp — 0.1 (a), 0.2 (b), 0.3 (c), and 0.4 (d). Right: adiabatic index for nuclear matter, including 
baryon, electron, positron and photon, with constant entropy S = 1, and proton fraction Yp = 0.1 (a), 0.2 (b), 0.3 (c), and 0.4 
(d). 



11 



10" 



Y p =0.2 
Y p =0.3 
Y p =0.4 
Y p =0.5 



(a) T=l MeV 



10 



10 



10 



10 

n [fm"] 



10" 



10 



to- 



10 



10 



10 



10 



10 



Y p =0.2 



Y p =0.4 



(b) T = 3.16 MeV 



10 



10" 



-5 -4 

10 10 
n B [fm ] 



10" 



10 



10 



____ Y p =0.1 








____ Y p =0.1 


i i 






Yp-0.2 








Yp-0.2 




/' 




Y p =0.3 








Y p =0.3 




i 
i 
1 




....... Yp-0.4 

- - - Yp-0.5 


1: ' ' '' 

li' w 

It ;! 




io 2 


....... Yp-0.4 

r Y p =0.5 




1 
i 

1 

1 

I 

1 




(c) T = 6.31 MeV 


Itf ;' 

1 / 

£•? J 






(d) T=10 MeV 




1 

1 
1 

i 

1 

I 


























. „1 




,i 







10"" 



10 

n B [fm""] 



10" 1 



10"" 



10 
n B [fm ' 



10"' 



FIG. 9: (Color online) Average mass number of heavy nuclei at temperatures of T = 1 (a), 3.16 (b), 6.31 (c), and 10 (d) MeV. 
The proton fraction ranges from Yp = 0.1 to 0.5. 
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FIG. 10: (Color online) Average proton number of heavy nuclei at temperatures of T = 1 (a), 3.16 (b), 6.31 (c), and 10 (d) 
MeV. The proton fraction ranges from Yp = 0.1 to 0.5. 



J. Comparisons with Lattimer-Swesty's EOS and 
H. Shen et aZ.'s EOS 

It is instructive to compare our EOS with the Lattimer- 
Swesty equation of state, that uses a simple liquid drop 
model with a non-relativistic Skyrme interaction, and the 
H. Shen et adequation of state, that uses the Thomas 
Fermi and variational approximations to a relativistic 
mean field model. The L-S EOS quoted in this work 
corresponds to the one with incompressibility coefficient 
A"=180 MeV. In Fig. Ho] the pressure for matter at T = 1, 
or 6.31 MeV and Y P = 0.05 or 0.4 is shown for our EOS, 
Lattimer-Swesty's (L-S) and H. Shen et al.'s (S-S) EOSs. 
Here the pressure includes contributions from electrons, 
positrons and photons. Below nuclear saturation density 
~ 0.16 fm -3 , ours and S-S EOS agree very well. The L- 
S EOS gives a slightly lower pressure at densities above 
10~ 3 fm~ 3 when T = 1 MeV and Y P = 0.05. Above 
saturation density, our EOS gives the largest pressure. 

Fig . 1 1 7l compares the mass fraction of heavy nuclei from 
our EOS, L-S EOS and S-S EOS, for matter at T = 1 or 



6.31 MeV, and Y P = 0.05 or 0.4. Although the pressure 
agrees well among the three EOSs, the mass fraction of 
heavy nuclei and alpha particles can be different and this 
contributes to the differences in entropy shown in Fig. 1181 
In addition our EOS may have a higher entropy in Fig llSI 
(b) because we may have smaller nuclei for n less than 
10~ 3 feu -3 compared to S-S and because we can have a 
distribution of several heavy nuclei, while the L-S and 
S-S EOSs use only a single average species. The aver- 
age A and Z of heavy nuclei is shown in Fig. Q1J] This 
figure clearly shows the effects of shell structure, that is 
included in our EOS but is neglected in both L-S and S-S 
EOSs. As a result, A and Z for our EOS have a serious 
of steps while A and Z for L-S and S-S EOSs increase 
smoothly with density. 

For T = 1 MeV and Y P = 0.05, the three EOSs agree 
very well on entropy, except that at low densities the L-S 
EOS gives a slightly smaller value. It is comprehensi- 
ble since the matter is dominated by free neutrons until 
very high density when heavy nuclei dominate (upper left 
panel of Fig. [T7|) . But heavy nuclei have low entropy at 
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low temperatures so the difference in entropy between 
the EOSs is small. For T = 1 MeV and Y P = 0.4, the 
differences in entropy between the three EOSs are larger. 
This is partially due to the smaller scale used. It is also 
because the different EOS predict different mass fractions 
of heavy nuclei, whose formation decreases the entropy. 
The L-S EOS has 80% heavy nuclei at 10" 6 fm" 3 (up- 
per right panel in Fig. [T7]) compared to 20% in our EOS 
and the S-S EOS, therefore the L-S EOS has the low- 
est entropy. At the higher temperature T = 6.31 MeV, 
the three EOSs agree well on entropy for different values 



of Yp, because free nucleons dominate at most densities 
(see lower panels in Fig. [TT]) . 

The differences in abundances and average A and Z 
between EOSs may arise because of differences in sym- 
metry energy and approximations made. For example 
L-S uses a very simple liquid drop model while the S-S 
EOS uses variational and Thomas Fermi approximations 
that are avoided in our EOS. These differences could be 
important for neutrino interactions and should be stud- 
ied further. 
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FIG. 11: (Color online) Mass fractions of free neutrons, free protons, alpha particles and heavy nuclei versus baryon density at 
four different temperatures T = 1 (a), 3.16 (b), 6.31 (c), and 10 (d) MeV, and fixed proton fraction Y p = 0.05. 



K. Adabatic compression tests 

One important test of thermodynamic consistency for 
our EOS table is to check that entropy is conserved as 
matter undergoes adiabatic compression. This is very 
closely related to the conservation of energy in the first 
law of thermodynamics. Let us consider the following 
adiabatic compression procudurc. 

1. Start from an energy per baryon E(n(to), T(io)), 
a density n(to), and a temperature T(to)- Here t 
stands for a sequence of times starting at to . The 
pressure P is then given by the EOS. 



2. Now compress the matter by increasing the density 
by an amount dn so that n(t\) = n(to) + dn with 
h =to + l. 

3. Next ensure that the first law of thermodynamics 
is satisfied by using a backward difference equation 
and requiring 

E(n(tx), T(ti)) = E{n{t Q ),T(t )) + P(n(h), T{t{))dn 

(10) 

by solving for the new pressure P(n(ii), T(t\)) and 
temperature T{t\). 

4. Finally the updated entropy S(n(t\), T(ti)) follows 
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FIG. 12: (Color online) Mass fractions of free neutrons, free protons, alpha particles and heavy nuclei versus baryon density at 
four different temperatures T = 1 (a), 3.16 (b), 6.31 (c), and 10 (d) MeV, and fixed proton fraction Y p = 0.2. 
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FIG. 13: (Color online) Mass fractions of free neutrons, free protons, alpha particles and heavy nuclei versus baryon density at 
four different temperatures T = 1 (a), 3.16 (b), 6.31 (c), and 10 (d) MeV, and fixed proton fraction Y p = 0.5. 
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FIG. 14: Energy (a), pressure (b) and proton fraction (c) versus density for zero temperature nuclear matter in chemical 
equilibrium. 
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FIG. 15: Mass vs Radius (a) and mass vs central density (b) for neutron star from zero temperature beta equilibrium nuclear 
matter as shown in Fig. Q3] 




FIG. 16: (Color online) Pressure versus baryon density from our EOS, Lattimer-Swesty's (L-S) and H. Shen et al.'s (S-S) EOS, 
with T = 1 MeV, Y p = 0.05 (a), T = 1 MeV, Y p = 0.4 (b), T = 6.31 MeV, Y p = 0.05 (c), and T = 6.31 MeV, Y p = 0.4 (d). 
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FIG. 17: (Color online) mass fraction of heavy nuclei versus baryon density from our EOS, Lattimer-Swesty's and H. Shen et 
al.'s EOS, with T = 1 MeV, Y p = 0.05 (a), T = 1 MeV, Y p = 0.4 (b), T = 6.31 MeV, Y p = 0.05 (c), and T = 6.31 MeV, Y p = 
0.4 (d). 




FIG. 18: (Color online) Entropy per baryon versus baryon density from our EOS, Lattimer-Swesty's and H. Shen et al.'s EOS, 
with T = 1 MeV, Y p = 0.05 (a), T = 1 MeV, Y p = 0.4 (b), T = 6.31 MeV, Y p = 0.05 (c), and T = 6.31 MeV, Y p = 0.4 (d). 
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FIG. 19: (Color online) Average mass number A (black upper curves) and atomic number Z (red lower curves) of heavy nuclei 
for our EOS (solid) , Lattimer-Swesty's (dot-dashed) and H. Shen et al.'s EOS (dashed) for T = 1 MeV, Y p = 0.4. 



from our EOS at the new density n{t\) and tem- 
perature T(ti). 

The entropy should be conserved in this adiabatic com- 
pression so that S is independent of t. In Fig. [20l the 
temperature and entropy versus density during adiabatic 
compression are shown for nuclear matter with fixed pro- 
ton fraction Yp = 0.15. The initial density is 5.16738 x 
1CP 6 fm -3 and initial temperature is 0.5 MeV. The solid 
curve is test result for our EOS table with lepton and 
photon contributions included in addition to contribu- 



tions from baryons. The dashed curve is a test result 
for our EOS table with baryons alone. The variation in 
entropy for both cases is less than 1%. For a comparison 
the test result for H. Shen EOS is also shown (as dot- 
dashed curve) for the similar initial condition. Note that 
it is important to use an accurate interpolation scheme 
with the EOS table to ensure that the first law is satisfied 
so that entropy is conserved. The adiabatic compression 
test result for H. Shen EOS was obtained using the rou- 
tine developed in Ref. [22j . 



IV. SUMMARY AND OUTLOOK 

In this paper we generate a new complete equation 
of state for use in supernova and neutron star merger 
simulations. We use a density dependent RMF model 
for nuclear matter at intermediate and high density with 
a spherical Wigner-Seitz approximation for nonuniform 
matter, which incorporates nuclear shell effects [3]. For 
nuclear matter at low density Q , we use a Virial expan- 
sion for a nonideal gas consisting of neutrons, protons, 
alpha particles and thousands of heavy nuclei from the 
finite range droplet model (FRDM) mass table [l8[. We 
include second order virial corrections for light elements 
A < 4, nuclear partition functions for heavy nuclei, and 



Coulomb corrections between electrons and heavy nuclei. 
As the density decreases, the mean field results match 
smoothly to the Virial gas. As the density goes down, 
the Virial expansion reduces to nuclear statistical equi- 
librium and the Virial expansion is exact in the low den- 
sity limit. The computed EOS covers 180,000 grid points 
in the temperature range T = ~ 80 MeV, the density 
range tib = 10~ 8 ~ 1.6 fm~ 3 , and the proton fraction 
range Yp = ~ 0.56. 

We use a hybrid interpolation scheme to generate a full 
EOS table on a fine grid that is thermodynamic consis- 
tent. This ensures that the first law of thermodynamics 
is satisfied and that entropy is conserved during adiabatic 
compression. Our EOS is an improvement over the ex- 
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FIG. 20: (Color online) Temperature (a) and entropy (b) versus density in the adiabatic compression for nuclear matter with 
fixed proton fraction 0.15. 



isting Lattimer-Swesty [TT| and H. Shen et al. fl2l [l3j. 
equations of state because our EOS includes thousands 
of heavy nuclei and is exact in the low density limit. 

We discuss the thermodynamic properties of our EOS 
in detail. The free energy per baryon, pressure, entropy 
per baryon and chemical potentials for neutrons, protons, 
and electron neutrinos are shown for matter at high and 
low temperatures and for large and small proton frac- 
tions. We also perform adiabatic compression tests for 
our EOS table and show that it preserves entropy to bet- 
ter than 1%. 

Finally, we show some comparisons between our EOS 
and two existing EOS tables. At high density, our RMF 
model reduces to the normal NL3 set. It produces a stiffer 
EOS at high density than the S-S EOS, which is already 
stiffer than the L-S EOS. One can find this difference, for 
example, in Fig.[T5]lor the pressure at high density. As a 
result, neutron stars obtained from our EOS have a large 
maximum mass of around 2.77 solar mass and radius of 
13.3 km. At low density, our model is essentially the 
virial expansion of a nonideal gas, which includes many 
nuclei. Our model usually has a broad distribution of nu- 
clei as shown in our previous paper [f| . This spread in the 
distribution of nuclei makes the composition of nuclear 
matter in our EOS different from that in the L-S and 
S-S EOS tables. This composition difference also intro- 
duces some differences in the entropy for nuclear matter 
at low temperature and large proton fraction, see for ex- 
ample the upper right panel in Fig. [TS] for matter with 



T = 1 MeV and Yp = 0.4. More importantly, the com- 
position of nuclear matter, especially the wide spreads in 
distribution of nuclei, could influence the position of the 
neutrinosphere and the spectra of radiated neutrinos and 
antineutrinos. 

The complete equation of state table is not only use- 
ful for core collapse supernova simulations, but also for 
other astrophysical applications, including the long term 
evolution of proto-neutron stars, and neutron star merg- 
ers. Our procedure, for matching mean field and virial 
calculations to produce a thermodynamically consistent 
EOS, can be used for a variety of effective interactions. In 
this paper we have presented an EOS based on the NL3 
interaction that is relatively stiff at high densities. In 
future work we will present additional EOS tables based 
on softer interactions that have lower pressures at high 
densities. 

Finally, our full EOS tables, both with and without 
lepton and photon contributions, are available for down- 
load as described in Appendix [Vj 

We thank Lorenz Huedepohl, Thomas Janka, Andreas 
Marek, Evan O'Connor, and Christian Ott for impor- 
tant help running astrophysical simulations to debug our 
equation of state. This work was supported in part by 
DOE grant DE-FG02-87ER40365. 
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V. APPENDIX: FORMAT OF EOS TABLES. 

Here we describe the physical quantities provided in 
the tables for the equation of state, NL3eosl.03.dat 
and NL3eosbl.03.dat and where they can be down- 
loaded. One should download the gzip compressed files 
NL3eosl.03.dat.gz and or NL3eosbl.03.dat.gz (that are 



about 100 MB each) and use gunzip to decompress them. 
The grid structures of these tables are indicated in Table 
ITTland contain approximately 517 MB of data each. The 
tables, a sample FORTRAN computer program, and a 
readme file are available for download both at our web- 
site 



http : / / cecelia. physics . indiana.edu/gang_shen_eos/ 



and, from the Electronic Physics Auxiliary Publication 
Service - EPAPS web site, 23]. Please check our web 
site for any updated information regarding these EOS 
tables. 

The tables contain a number of quantities. For a single 
triplet of T, tib, and Yp, there are 16 items that make up 
a row of the table. These items are described below. In 
the table NL3eosbl.03.dat, only the contribution from 
baryons is taken into account for items 4,5,6. In the 
table NL3eosl.03.dat, the contributions from electrons, 
positrons, and photons are also included. The electron 
mass is 0.511 MeV. 



A. Items in EOS tables. 

1. Temperature T [MeV]. The range of temperature 
is first T=0 then from 10" 8 - 10 1 - 875 MeV, see 
Table III 

2. Proton fraction Yp. The range of proton fraction 
is first 0, and then from, 0.05 ~ 0.56. The step in 
proton fraction is 0.01, see Table [TT1 

3. Baryon number density n [fin -3 ]. The range of den- 
sity is from 10~ 8 to 10 0175 fm" 3 . See Table__] 

4. Free energy per baryon F/A [MeV] which has sub- 
tracted the free nucleon mass 939 MeV. 

5. Pressure P [MeV/fm 3 ]. 

6. Entropy per baryon S/A [kg]. 

7. Chemical potential for neutrons \x n [MeV]. The tab- 
ulated value is relative to the nucleon mass 939 
MeV. 

8. Chemical potential for protons fi p [MeV]. The tab- 
ulated value is relative to the nucleon mass. 



9. Chemical potential for electrons ^ e [MeV]. The tab- 
ulated value includes the electron mass. 

10. Average mass number A of heavy nuclei with A > 
4, which exclude alpha particles. 

11. Average proton number Z of heavy nuclei with A > 
4, which exclude alpha particles. 

12. Mass fraction of free neutrons. 

13. Mass fraction of free protons. 

14. Mass fraction of alpha particles. 

15. Mass fraction of heavy nuclei with A > 4, which 
exclude alpha particles. 

16. Effective nucleon mass M* [MeV]. In uniform mat- 
ter it is obtained from RMF theory. For virial gas 
and non-uniform matter, it is chosen to be the free 
nucleon mass M* = 939 MeV. 



B. Sample FORTRAN computer program readeos.f 

The Fortran program readeos.f includes a very short 
main program that calls the subroutine load_table, to 
read NL3eosl.03.dat or NL3eosbl.03.dat, and then calls 
the subroutine readeos with inputs T (in MeV), proton 
fraction Y p , and density n (in fm -3 ). The subroutine 
readeos uses triliner interpolation (in T, Y p , and n) to 
return the above 16 values plus the internal energy per 
baryon (in MeV) and the chemical potential for electron 
neutrinos in chemical equilibrium (in MeV). Note that 
one needs to call load_table only once and then one can 
call readeos many times. For further details please see 
the comments in readeos.f. 
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